BioGSP: Spectral Graph Wavelet Transform for Spatial Data (New Workflow)

Introduction

BioGSP provides Spectral Graph Wavelet Transform (SGWT) analysis for spatial biological data. This package enables multi-scale analysis of spatial patterns using graph signal processing techniques.

Key Capabilities

• Object-oriented SGWT workflow with clear separation of concerns
• Multi-scale spatial analysis using spectral graph wavelets
• Flexible data input with custom column naming
• Energy-normalized similarity analysis with comprehensive metrics
• Multiple wavelet kernels (Mexican Hat, Meyer, Heat)
• Comprehensive visualization tools
• Detailed inverse transform results (low-pass, band-pass approximations)

New Workflow Overview

The new workflow consists of four main steps:

1. initSGWT(): Initialize SGWT object with data and parameters
2. runSpecGraph(): Build spectral graph structure
3. runSGWT(): Perform forward and inverse SGWT transforms
4. runSGCC(): Calculate energy-normalized weighted similarity

Setup

# Load required packages
library(ggplot2)
library(patchwork)
library(viridis)
library(Matrix)
library(igraph)
library(RANN)
library(RSpectra)

# Load BioGSP functions
source('../R/simulation.R')
source('../R/sgwt_core.R')
source('../R/sgwt_main.R')
source('../R/utilities.R')
source('../R/visualization.R')

# Set random seed for reproducibility
set.seed(123)

New SGWT Workflow Demonstration

Generate Example Pattern

# Create a spatial pattern with concentric circles
demo_pattern <- simulate_multiscale(
  grid_size = 40,        # Moderate size for fast computation
  n_centers = 1,         # Single center pattern
  Ra_seq = 5,            # Inner circle radius
  Rb_seq = 10,           # Outer ring radius
  seed = 123
)[[1]]

## Generating 1 multiscale patterns...
## 
  |                                                                            
  |                                                                      |   0%
  |                                                                            
  |======================================================================| 100%
## 
## Multiscale simulation completed!

# Display pattern structure
cat("Generated pattern dimensions:", nrow(demo_pattern), "x", ncol(demo_pattern), "\n")

## Generated pattern dimensions: 1600 x 4

cat("Column names:", paste(colnames(demo_pattern), collapse = ", "), "\n")

## Column names: X, Y, signal_1, signal_2

cat("Signals available:", sum(demo_pattern$signal_1), "inner pixels,", sum(demo_pattern$signal_2), "outer pixels\n")

## Signals available: 81 inner pixels, 628 outer pixels

Visualize Input Pattern

# Create visualizations for both signals
p1 <- ggplot(demo_pattern, aes(X, Y, fill = factor(signal_1))) +
  geom_tile() +
  scale_fill_manual(values = c("0" = "white", "1" = "#e31a1c"), 
                    name = "Signal", labels = c("Background", "Inner Circle")) +
  coord_fixed() + theme_void() + ggtitle("Signal 1: Inner Circles")

p2 <- ggplot(demo_pattern, aes(X, Y, fill = factor(signal_2))) +
  geom_tile() +
  scale_fill_manual(values = c("0" = "white", "1" = "#1f78b4"), 
                    name = "Signal", labels = c("Background", "Outer Ring")) +
  coord_fixed() + theme_void() + ggtitle("Signal 2: Outer Rings")

print(p1 + p2)

Step 1: Initialize SGWT Object

# Initialize SGWT object with custom column names
SG <- initSGWT(
  data.in = demo_pattern,
  x_col = "X",           # Custom X coordinate column name
  y_col = "Y",           # Custom Y coordinate column name
  signals = c("signal_1", "signal_2"),  # Analyze both signals
  J = 4,                 # Number of wavelet scales
  scaling_factor = 2,    # Scale progression factor
  kernel_type = "heat"
)

# Display initialized object
print(SG)

## SGWT Object
## ===========
## Data:
##   Dimensions: 1600 x 4 
##   Coordinates: X , Y 
##   Signals: signal_1, signal_2 
## 
## Parameters:
##   k (neighbors): heat 
##   J (scales): 4 
##   Kernel type: heat 
##   Laplacian type: 
## 
## Status:
##   Graph computed: FALSE 
##   Forward computed: FALSE 
##   Inverse computed: FALSE

Step 2: Build Spectral Graph

# Build spectral graph structure
SG <- runSpecGraph(SG, k = 12, laplacian_type = "normalized", length_eigenvalue = 100, verbose = TRUE)

## Building graph from spatial coordinates...
## Computing Laplacian and eigendecomposition...
## Auto-generated scales: 0.3761, 0.1881, 0.094, 0.047 
## Graph construction completed.

# Check updated object
cat("Graph construction completed!\n")

## Graph construction completed!

cat("Adjacency matrix dimensions:", dim(SG$Graph$adjacency_matrix), "\n")

## Adjacency matrix dimensions: 1600 1600

cat("Number of eigenvalues computed:", length(SG$Graph$eigenvalues), "\n")

## Number of eigenvalues computed: 100

Visualize Fourier Modes (Eigenvectors)

# Visualize Fourier modes (eigenvectors) - 5 low and 5 high frequency modes
fourier_modes <- plot_FM(SG, mode_type = "both", n_modes = 5, ncol = 5)

cat("Fourier Modes Visualization:\n")

## Fourier Modes Visualization:

cat("- Low-frequency modes (2-6): Smooth spatial patterns, excluding DC component\n")

## - Low-frequency modes (2-6): Smooth spatial patterns, excluding DC component

cat("- High-frequency modes (96-100): Fine-detailed spatial patterns\n")

## - High-frequency modes (96-100): Fine-detailed spatial patterns

cat("- Each mode shows its eigenvalue (λ) indicating frequency content\n")

## - Each mode shows its eigenvalue (λ) indicating frequency content

Step 3: Run SGWT Analysis

# Perform SGWT forward and inverse transforms
SG <- runSGWT(SG, verbose = TRUE)

## Performing SGWT analysis for 2 signals...
## Using batch processing for efficiency...
## SGWT analysis completed.

# Display final object with all results
print(SG)

## SGWT Object
## ===========
## Data:
##   Dimensions: 1600 x 4 
##   Coordinates: X , Y 
##   Signals: signal_1, signal_2 
## 
## Parameters:
##   k (neighbors): heat 
##   J (scales): 4 
##   Kernel type: heat 
##   Laplacian type: 
##   Scales: 0.3761, 0.1881, 0.094, 0.047 
## 
## Status:
##   Graph computed: TRUE 
##   Forward computed: TRUE 
##   Inverse computed: TRUE 
## 
## Reconstruction Errors:
##    signal_1 : 0.175508 
##    signal_2 : 0.375262

Visualize SGWT Decomposition

# Plot SGWT decomposition results for signal_1
decomp_plots <- plot_sgwt_decomposition(
  SG = SG,
  signal_name = "signal_1",
  plot_scales = 1:3,     # Show first 3 scales
  ncol = 3
)

Energy Analysis

# Analyze energy distribution across scales for signal_1
energy_analysis <- sgwt_energy_analysis(SG, "signal_1")
print(energy_analysis)

##       scale    energy energy_ratio scale_value   signal
## 1  low_pass 46.351714   0.75046488  0.37613716 signal_1
## 2 wavelet_1  1.647846   0.02667971  0.37613716 signal_1
## 3 wavelet_2  3.510442   0.05683638  0.18806858 signal_1
## 4 wavelet_3  5.180645   0.08387806  0.09403429 signal_1
## 5 wavelet_4  5.073355   0.08214096  0.04701714 signal_1

# Create energy distribution plot
energy_plot <- ggplot(energy_analysis, aes(x = scale, y = energy_ratio)) +
  geom_bar(stat = "identity", fill = "steelblue", alpha = 0.7) +
  geom_text(aes(label = paste0(round(energy_ratio * 100, 1), "%")), 
            vjust = -0.5, size = 3) +
  labs(title = "Energy Distribution Across SGWT Scales (Signal 1)",
       x = "Scale", y = "Energy Ratio") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

print(energy_plot)

Advanced Features Demonstration

Different Wavelet Kernels

# Compare different kernel types
kernels <- c("mexican_hat", "meyer", "heat")
kernel_results <- list()

for (kn in kernels) {
  # Initialize with different kernel
  SG_temp <- initSGWT(
    data.in = demo_pattern,
    x_col = "X", y_col = "Y",
    signals = "signal_1",
    J = 3,                    # Reduced for faster computation
    kernel_type = kn
  )
  
  # Run full workflow
  SG_temp <- runSpecGraph(SG_temp, k = 12, laplacian_type = "normalized", length_eigenvalue = 100,verbose = FALSE)
  SG_temp <- runSGWT(SG_temp, verbose = FALSE)
  
  kernel_results[[kn]] <- SG_temp
}

# Compare reconstruction errors
comparison_df <- data.frame(
  Kernel = kernels,
  RMSE = sapply(kernel_results, function(x) x$Inverse$signal_1$reconstruction_error),
  stringsAsFactors = FALSE
)

print("Kernel Performance Comparison:")

## [1] "Kernel Performance Comparison:"

print(comparison_df)

##                  Kernel      RMSE
## mexican_hat mexican_hat 0.1579043
## meyer             meyer 0.1720928
## heat               heat 0.1329538

# Plot comparison
ggplot(comparison_df, aes(x = Kernel, y = RMSE, fill = Kernel)) +
  geom_bar(stat = "identity", alpha = 0.7) +
  geom_text(aes(label = round(RMSE, 6)), vjust = -0.5) +
  scale_fill_viridis_d() +
  labs(title = "SGWT Reconstruction Error by Kernel Type",
       x = "Kernel Type", y = "RMSE") +
  theme_minimal()

Step 4: Pattern Similarity Analysis

# Calculate similarity between signal_1 and signal_2 in same object
similarity_within <- runSGCC("signal_1", "signal_2", SG = SG, return_parts = TRUE)

cat("Pattern Similarity Analysis (within same graph):\n")

## Pattern Similarity Analysis (within same graph):

cat(sprintf("Overall similarity: %.4f\n", similarity_within$S))

## Overall similarity: -0.1200

cat(sprintf("Low-frequency similarity: %.4f\n", similarity_within$c_low))

## Low-frequency similarity: 0.0005

cat(sprintf("Non-low-frequency similarity: %.4f\n", similarity_within$c_nonlow))

## Non-low-frequency similarity: -0.5266

cat(sprintf("Energy weights - Low: %.3f, Non-low: %.3f\n", 
           similarity_within$w_low, similarity_within$w_NL))

## Energy weights - Low: 0.771, Non-low: 0.229

# Generate a second pattern for cross-comparison
pattern_2 <- simulate_multiscale(
  grid_size = 40,
  n_centers = 1,
  Ra_seq = 7,            # Different inner radius
  Rb_seq = 12,           # Different outer radius
  seed = 456
)[[1]]

## Generating 1 multiscale patterns...
## 
  |                                                                            
  |                                                                      |   0%
  |                                                                            
  |======================================================================| 100%
## 
## Multiscale simulation completed!

# Create second SGWT object
SG2 <- initSGWT(pattern_2, x_col = "X", y_col = "Y", signals = "signal_1", 
                J = 4)
SG2 <- runSpecGraph(SG2, k = 12, laplacian_type = "normalized",length_eigenvalue = 100, verbose = FALSE)
SG2 <- runSGWT(SG2, verbose = FALSE)

# Calculate cross-object similarity
similarity_cross <- runSGCC(SG, SG2, return_parts = TRUE)

cat("\nCross-object Pattern Similarity:\n")

## 
## Cross-object Pattern Similarity:

cat(sprintf("Overall similarity: %.4f\n", similarity_cross$S))

## Overall similarity: 0.9153

cat(sprintf("Low-frequency similarity: %.4f\n", similarity_cross$c_low))

## Low-frequency similarity: 0.9355

cat(sprintf("Non-low-frequency similarity: %.4f\n", similarity_cross$c_nonlow))

## Non-low-frequency similarity: 0.8459

# Visualize both patterns for comparison
p_comp1 <- ggplot(demo_pattern, aes(X, Y, fill = factor(signal_1))) +
  geom_tile() + scale_fill_manual(values = c("0" = "white", "1" = "#e31a1c")) +
  coord_fixed() + theme_void() + theme(legend.position = "none") +
  ggtitle("Pattern A (R_in=5, R_out=10)")

p_comp2 <- ggplot(pattern_2, aes(X, Y, fill = factor(signal_1))) +
  geom_tile() + scale_fill_manual(values = c("0" = "white", "1" = "#1f78b4")) +
  coord_fixed() + theme_void() + theme(legend.position = "none") +
  ggtitle("Pattern B (R_in=7, R_out=12)")

print(p_comp1 + p_comp2)

Low-frequency Only Analysis

# Compare using only low-frequency components
similarity_low <- runSGCC("signal_1", "signal_2", SG = SG, low_only = TRUE, return_parts = TRUE)

cat("Low-frequency Only Similarity Analysis:\n")

## Low-frequency Only Similarity Analysis:

cat(sprintf("Low-frequency similarity: %.4f\n", similarity_low$c_low))

## Low-frequency similarity: 0.0005

cat(sprintf("Overall score (same as low-freq): %.4f\n", similarity_low$S))

## Overall score (same as low-freq): 0.0005

cat("Note: Non-low-frequency components are ignored (c_nonlow = NA)\n")

## Note: Non-low-frequency components are ignored (c_nonlow = NA)

# Compare with full analysis
cat("\nComparison:\n")

## 
## Comparison:

cat(sprintf("Full analysis similarity: %.4f\n", similarity_within$S))

## Full analysis similarity: -0.1200

cat(sprintf("Low-only similarity: %.4f\n", similarity_low$S))

## Low-only similarity: 0.0005

cat(sprintf("Difference: %.4f\n", abs(similarity_within$S - similarity_low$S)))

## Difference: 0.1205

Multiple Signal Analysis

# Analyze energy distribution for both signals
energy_signal1 <- sgwt_energy_analysis(SG, "signal_1")
energy_signal2 <- sgwt_energy_analysis(SG, "signal_2")

# Combine for comparison
energy_comparison <- rbind(energy_signal1, energy_signal2)

# Plot comparison
ggplot(energy_comparison, aes(x = scale, y = energy_ratio, fill = signal)) +
  geom_bar(stat = "identity", position = "dodge", alpha = 0.7) +
  geom_text(aes(label = paste0(round(energy_ratio * 100, 1), "%")), 
            position = position_dodge(width = 0.9), vjust = -0.5, size = 3) +
  scale_fill_viridis_d() +
  labs(title = "Energy Distribution Comparison: Signal 1 vs Signal 2",
       x = "Scale", y = "Energy Ratio", fill = "Signal") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Similarity Space Visualization with Multiple Patterns

# Generate 9 paired patterns using simulate_multiscale
cat("Generating 9 paired patterns for similarity analysis...\n")

## Generating 9 paired patterns for similarity analysis...

patterns_9 <- simulate_multiscale(
  grid_size = 40,
  n_centers = 1,
  Ra_seq = c(1, 5, 10),    # Inner circle radii
  Rb_seq = c(2, 6, 15),    # Outer ring radii
  seed = 123
)

## Generating 9 multiscale patterns...
## 
  |                                                                            
  |                                                                      |   0%
  |                                                                            
  |========                                                              |  11%
  |                                                                            
  |================                                                      |  22%
  |                                                                            
  |=======================                                               |  33%
  |                                                                            
  |===============================                                       |  44%
  |                                                                            
  |=======================================                               |  56%
  |                                                                            
  |===============================================                       |  67%
  |                                                                            
  |======================================================                |  78%
  |                                                                            
  |==============================================================        |  89%
  |                                                                            
  |======================================================================| 100%
## 
## Multiscale simulation completed!

cat("Generated", length(patterns_9), "patterns\n")

## Generated 9 patterns

cat("Pattern names:", names(patterns_9), "\n")

## Pattern names: simulated_Ra_1_Rb_2 simulated_Ra_1_Rb_6 simulated_Ra_1_Rb_15 simulated_Ra_5_Rb_2 simulated_Ra_5_Rb_6 simulated_Ra_5_Rb_15 simulated_Ra_10_Rb_2 simulated_Ra_10_Rb_6 simulated_Ra_10_Rb_15

# Visualize all 9 patterns first
pattern_viz <- visualize_multiscale(
  patterns_9, 
  Ra_seq = c(1, 5, 10), 
  Rb_seq = c(2, 6, 15),
  show_title = TRUE
)
print(pattern_viz)

# Create SGWT objects for all 9 patterns and compute similarities
similarity_results <- list()
sgwt_objects <- list()

cat("\nProcessing patterns and computing SGWT analysis...\n")

## 
## Processing patterns and computing SGWT analysis...

# Process each pattern
for (i in seq_along(patterns_9)) {
  pattern_name <- names(patterns_9)[i]
  pattern_data <- patterns_9[[i]]
  
  cat("Processing", pattern_name, "...\n")
  
  # Create SGWT object
  SG_temp <- initSGWT(
    data.in = pattern_data,
    x_col = "X", 
    y_col = "Y", 
    signals = c("signal_1", "signal_2"),
    J = 4,
    kernel_type = "heat"
  )
  
  # Build graph and run SGWT
  SG_temp <- runSpecGraph(SG_temp, k = 12, laplacian_type = "normalized", 
                         length_eigenvalue = 50, verbose = FALSE)
  SG_temp <- runSGWT(SG_temp, verbose = FALSE)
  
  sgwt_objects[[pattern_name]] <- SG_temp
  
  # Compute within-pattern similarity (signal_1 vs signal_2)
  sim_within <- runSGCC("signal_1", "signal_2", SG = SG_temp, return_parts = TRUE)
  similarity_results[[pattern_name]] <- sim_within
}

## Processing simulated_Ra_1_Rb_2 ...
## Processing simulated_Ra_1_Rb_6 ...
## Processing simulated_Ra_1_Rb_15 ...
## Processing simulated_Ra_5_Rb_2 ...
## Processing simulated_Ra_5_Rb_6 ...
## Processing simulated_Ra_5_Rb_15 ...
## Processing simulated_Ra_10_Rb_2 ...
## Processing simulated_Ra_10_Rb_6 ...
## Processing simulated_Ra_10_Rb_15 ...

# Print similarity results summary
cat("\nSimilarity Results Summary (signal_1 vs signal_2 within each pattern):\n")

## 
## Similarity Results Summary (signal_1 vs signal_2 within each pattern):

similarity_df <- data.frame(
  Pattern = names(similarity_results),
  S = sapply(similarity_results, function(x) x$S),
  c_low = sapply(similarity_results, function(x) x$c_low),
  c_nonlow = sapply(similarity_results, function(x) x$c_nonlow),
  w_low = sapply(similarity_results, function(x) x$w_low),
  w_NL = sapply(similarity_results, function(x) x$w_NL),
  stringsAsFactors = FALSE
)

print(similarity_df)

##                                     Pattern          S       c_low   c_nonlow
## simulated_Ra_1_Rb_2     simulated_Ra_1_Rb_2  0.9938194  0.99422031  0.9931752
## simulated_Ra_1_Rb_6     simulated_Ra_1_Rb_6  0.8273410  0.85675482  0.7693110
## simulated_Ra_1_Rb_15   simulated_Ra_1_Rb_15  0.2081565  0.28593433  0.0382921
## simulated_Ra_5_Rb_2     simulated_Ra_5_Rb_2  0.6962827  0.78080645  0.5001333
## simulated_Ra_5_Rb_6     simulated_Ra_5_Rb_6  0.2252059  0.38358293 -0.1237262
## simulated_Ra_5_Rb_15   simulated_Ra_5_Rb_15 -0.3014023 -0.21588448 -0.4877800
## simulated_Ra_10_Rb_2   simulated_Ra_10_Rb_2  0.2738171  0.36207226  0.0608327
## simulated_Ra_10_Rb_6   simulated_Ra_10_Rb_6 -0.1695127 -0.07201739 -0.3860797
## simulated_Ra_10_Rb_15 simulated_Ra_10_Rb_15 -0.9585733 -0.96326437 -0.9444922
##                           w_low      w_NL
## simulated_Ra_1_Rb_2   0.6163359 0.3836641
## simulated_Ra_1_Rb_6   0.6636265 0.3363735
## simulated_Ra_1_Rb_15  0.6859265 0.3140735
## simulated_Ra_5_Rb_2   0.6988533 0.3011467
## simulated_Ra_5_Rb_6   0.6878095 0.3121905
## simulated_Ra_5_Rb_15  0.6854755 0.3145245
## simulated_Ra_10_Rb_2  0.7070267 0.2929733
## simulated_Ra_10_Rb_6  0.6895670 0.3104330
## simulated_Ra_10_Rb_15 0.7501058 0.2498942

# Extract Ra and Rb values for better labeling
similarity_df$Ra <- as.numeric(gsub(".*Ra_([0-9.]+)_Rb.*", "\\1", similarity_df$Pattern))
similarity_df$Rb <- as.numeric(gsub(".*Rb_([0-9.]+)$", "\\1", similarity_df$Pattern))
similarity_df$Label <- paste0("Ra=", similarity_df$Ra, ",Rb=", similarity_df$Rb)

# Create similarity space visualization
similarity_plot <- visualize_similarity_xy(
  similarity_results,
  point_size = 4,
  point_color = "steelblue",
  add_diagonal = TRUE,
  add_axes_lines = TRUE,
  title = "SGWT Similarity Space: 9 Paired Patterns (signal_1 vs signal_2)",
  show_labels = TRUE
)

print(similarity_plot)

# Create enhanced plot with color coding by Ra values
library(ggplot2)
plot_data <- data.frame(
  c_low = similarity_df$c_low,
  c_nonlow = similarity_df$c_nonlow,
  Ra = factor(similarity_df$Ra),
  Rb = factor(similarity_df$Rb),
  Label = similarity_df$Label,
  S = similarity_df$S
)

enhanced_plot <- ggplot(plot_data, aes(x = c_low, y = c_nonlow, color = Ra, shape = Rb)) +
  geom_point(size = 5, alpha = 0.8) +
  geom_text(aes(label = Label), vjust = -0.8, hjust = 0.5, size = 3, color = "black") +
  xlim(-1, 1) + ylim(-1, 1) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "gray60", alpha = 0.7) +
  geom_vline(xintercept = 0, linetype = "dashed", color = "gray60", alpha = 0.7) +
  geom_abline(slope = 1, intercept = 0, linetype = "dotted", color = "gray40", alpha = 0.7) +
  geom_abline(slope = -1, intercept = 0, linetype = "dotted", color = "gray40", alpha = 0.7) +
  scale_color_viridis_d(name = "Inner Radius (Ra)") +
  scale_shape_manual(values = c(16, 17, 15), name = "Outer Radius (Rb)") +
  labs(
    title = "Enhanced Similarity Space: Color by Ra, Shape by Rb",
    x = "Low-frequency Similarity (c_low)",
    y = "Non-low-frequency Similarity (c_nonlow)"
  ) +
  theme_minimal() +
  theme(
    plot.title = element_text(hjust = 0.5, size = 14, face = "bold"),
    legend.position = "right"
  )

print(enhanced_plot)

# Analysis of patterns
cat("\nPattern Analysis:\n")

## 
## Pattern Analysis:

cat("Ra (Inner Radius) Effect:\n")

## Ra (Inner Radius) Effect:

for (ra in unique(similarity_df$Ra)) {
  subset_data <- similarity_df[similarity_df$Ra == ra, ]
  cat(sprintf("  Ra=%g: Mean c_low=%.3f, Mean c_nonlow=%.3f, Mean S=%.3f\n", 
              ra, mean(subset_data$c_low), mean(subset_data$c_nonlow), mean(subset_data$S)))
}

##   Ra=1: Mean c_low=0.712, Mean c_nonlow=0.600, Mean S=0.676
##   Ra=5: Mean c_low=0.316, Mean c_nonlow=-0.037, Mean S=0.207
##   Ra=10: Mean c_low=-0.224, Mean c_nonlow=-0.423, Mean S=-0.285

cat("\nRb (Outer Radius) Effect:\n")

## 
## Rb (Outer Radius) Effect:

for (rb in unique(similarity_df$Rb)) {
  subset_data <- similarity_df[similarity_df$Rb == rb, ]
  cat(sprintf("  Rb=%g: Mean c_low=%.3f, Mean c_nonlow=%.3f, Mean S=%.3f\n", 
              rb, mean(subset_data$c_low), mean(subset_data$c_nonlow), mean(subset_data$S)))
}

##   Rb=2: Mean c_low=0.712, Mean c_nonlow=0.518, Mean S=0.655
##   Rb=6: Mean c_low=0.389, Mean c_nonlow=0.087, Mean S=0.294
##   Rb=15: Mean c_low=-0.298, Mean c_nonlow=-0.465, Mean S=-0.351

cat("\nSimilarity Space Interpretation:\n")

## 
## Similarity Space Interpretation:

cat("- Each point represents similarity between signal_1 (inner circle) and signal_2 (outer ring)\n")

## - Each point represents similarity between signal_1 (inner circle) and signal_2 (outer ring)

cat("- Color indicates inner radius (Ra): darker = smaller radius\n")

## - Color indicates inner radius (Ra): darker = smaller radius

cat("- Shape indicates outer radius (Rb): circle/triangle/square for different sizes\n")

## - Shape indicates outer radius (Rb): circle/triangle/square for different sizes

cat("- Position shows frequency-domain similarity characteristics\n")

## - Position shows frequency-domain similarity characteristics

Advanced Workflow Features

Custom Parameters

# Demonstrate advanced parameter customization
SG_advanced <- initSGWT(
  data.in = demo_pattern,
  x_col = "X", y_col = "Y",
  signals = "signal_1",
  J = 6,                     # More scales for finer analysis
  scaling_factor = 1.5,      # Closer scales
  kernel_type = "heat"       # Heat kernel
)

SG_advanced <- runSpecGraph(SG_advanced, k = 15, laplacian_type = "randomwalk",length_eigenvalue = 30, verbose = FALSE)
SG_advanced <- runSGWT(SG_advanced, verbose = FALSE)

cat("Advanced Parameters Results:\n")

## Advanced Parameters Results:

cat("Number of scales:", length(SG_advanced$Parameters$scales), "\n")

## Number of scales: 6

cat("Scales:", paste(round(SG_advanced$Parameters$scales, 4), collapse = ", "), "\n")

## Scales: 0.1226, 0.0817, 0.0545, 0.0363, 0.0242, 0.0161

cat("Reconstruction error:", round(SG_advanced$Inverse$signal_1$reconstruction_error, 6), "\n")

## Reconstruction error: 0.256808

Summary

Key Takeaways

1. New Object-Oriented Workflow: Clear separation of initialization, graph building, analysis, and similarity computation
2. Flexible Input: Supports custom column naming for coordinates and signals
   
3. Multiple Kernels: Choose from Mexican Hat, Meyer, or Heat kernel families
4. Comprehensive Similarity: Energy-normalized weighted similarity with detailed diagnostics
5. Advanced Analysis: Low-frequency only analysis, cross-object comparison, multiple signal support

New Workflow Benefits

# Complete workflow in 4 clear steps:
SG <- initSGWT(data, signals = c("signal1", "signal2"))  # 1. Initialize
SG <- runSpecGraph(SG, k = 12, laplacian_type = "normalized")  # 2. Build graph  
SG <- runSGWT(SG)                                        # 3. Run SGWT
similarity <- runSGCC("signal1", "signal2", SG = SG)    # 4. Compare patterns

Parameter Guidelines

• Small datasets (<200 points): k=8-12, J=3-4
• Medium datasets (200-1000 points): k=12-20, J=4-6
  
• Large datasets (>1000 points): k=15-25, J=5-8

Next Steps

• Explore your own spatial datasets using the demonstrated workflow
• Experiment with different kernel types and parameters
• Use comprehensive similarity analysis for comparative studies
• Refer to function documentation for detailed parameter specifications

---

BioGSP Package - Biological Graph Signal Processing for Spatial Data Analysis